Chaotic dynamics in two-dimensional Rayleigh-Benard convection 
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We investigate the origin of various convective patterns using bifurcation diagrams that are con- 
structed using direct numerical simulations. We perform two-dimensional pseudospectral simula- 
tions for a Prandtl number 6.8 fluid that is confined in a box with aspect ratio T = 2y2. Steady 
convective rolls are born from the conduction state through a pitchfork bifurcation at r = 1, where 
r is the reduced Rayleigh number. These fixed points bifurcate successively to time-periodic and 
quasiperiodic rolls through Hopf and Neimark-Sacker bifurcations at r ~ 80 and r ~ 500 respec- 
^" tively. The system becomes chaotic at r ~ 750 through a quasiperiodic route to chaos. The size of 

the chaotic attractor increases at r ~ 840 through an "attractor-merging crisis" which also results 
in travelling chaotic rolls. We also observe coexistence of stable fixed points and a chaotic attractor 
for 846 < r < 849 as a result of a subcritical Hopf bifurcation. Subsequently the chaotic attrac- 
tor disappears through a "boundary crisis" and only stable fixed points remain. Later these fixed 
rf\ . points become periodic and chaotic through another set of bifurcations which ultimately leads to 

turbulence. 

Ch ' PACS numbers: 47.20.Bp, 47.27.ek, 47.52. +j 
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4j ! I. INTRODUCTION 

a: 

jy! . Rayleigh-Benard convection (RBC), an idealized version of thermal convection, is studied for understanding pattern- 

forming instabilities, chaos, and turbulence in atmosphere, astrophysics, crystal growth, etc. |l|. The flow dynamics 
in RBC is governed by two nondimensional parameters: Rayleigh number R, which is the ratio of buoyancy and the 

^*y dissipative terms, and the Prandtl number P, which is the ratio of the kinematic viscosity and the thermal diffusivity. 
rS^ Convection starts at the critical Rayleigh number R c that has been found to be independent of the Prandtl number. 

Oh The numerical value of R c is approximately 1708 for no-slip boundary conditions, and 657.5 for free-slip boundary 

conditions [l]. 
^-H ■ Secondary instabilities and the resulting patterns of RBC depend quite critically on the Prandtl number. For low 

£> ' Prandtl-number (low-P) and zero Prandtl number (zero-P) convection, the inertial term u • Vu is quite important 
and it generates vertical vorticity. Consequently the flow pattern becomes three-dimensional, and oscillatory waves 
along the horizontal axes are generated just near the onset of convection. For large Prandtl number (large-P) con- 
vection, however, vertical vorticity is absent near the onset and two-dimensional (2D) rolls survive till large Rayleigh 
numbers [lj, 0] . This feature of large-P convection makes 2D simulations quite relevant for this regime. 

Krishnamurti [3[ carried out extensive convection experiments on fluids of different Prandtl numbers (0.025 < P < 
8500) and studied various convective states such as two-dimensional, time-periodic, and chaotic rolls, and turbulent 
structures. Gollub and Benson [4] performed RBC experiments on water for P — 2.5 and 5.0 (varied by changing the 
mean temperature) and aspect ratios T — 2.4 and 3.5. They observed a variety of convective states and several routes 
to chaos. For T = 3.5 and P = 5.0, they found a quasiperiodic state at the reduced Rayleigh number r — R/R c — 32 
and subsequent phase locking at r ~ 44.4. The flow became chaotic at around r ~ 46, thus exhibiting a quasiperiodic 
route to chaos. Gollub and Benson [4{ also observed a period-doubling route to chaos for T — 3.5 and P = 2.5, and 
a quasiperiodic route to chaos with three frequencies (Ruelle-Takens) for T = 2.4 and P — 5.0. Intermittent chaos 
too was observed in their experiments. Maurer and Libchaber [5J observed quasiperiodic and frequency- locking route 
to chaos in helium (P ~ 0.7). Experiments of Giglio et al. [6j on the convection of water (P ~ 7) showed a period 
doubling route to chaos. Berge et al. Q observed intermittency in their RBC experiments on silicon oil (P ~ 100). 
Libchaber et al. [8[ studied convection in mercury (P ~ 0.02) in the presence of an applied mean magnetic field and 
observed generation of chaos through period doubling and quasiperiodic routes. 

Direct numerical simulations (DNS) are used extensively to study convection. Curry et al. [9( performed DNS in 
three dimensions for a P = 10 fluid under free-slip boundary conditions. They reported steady convection till r ~ 40, 
after which a limit cycle is observed till r ~ 45. Subsequently they reported quasiperiodicity (r ~ 45 — 55), phase 
locking (r ~ 55 — 65), and chaos (r > 65). Yahata [1(| observed a similar route to chaos in his DNS for P = 5 
, T x — 3.5 and T y = 2. Mukutmoni and Yang [llj performed DNS for convection in a rectangular enclosure with 
insulated side walls and found a period-two state, but a quasiperiodic route to chaos. Yahata [l(| and Nishikawa 
and Yahata [l2| performed direct numerical simulations of RBC for the setup similar to Gollub and Benson's [J] 
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experiments and observed quasiperiodic route to chaos. Yahata [10| also observed frequency locking. Gelfgat [13[ 
numerically studied the effects of the aspect ratio on the critical Rayleigh number and the instability modes. 

Three-dimensional DNS are quite expensive computationally, so a large number of two-dimensional DNS have been 
performed. Moore and Weiss 14] simulated for P = 6.8 using a finite-difference method under free-slip boundary 
conditions and reported that the Nusselt number Nu sw 2r 1 l' i for 5 < r < P 15 , and Nu sw 2r 0365 for higher Rayleigh 
numbers. McLaughlin and Orszag [15| simulated RBC in air (P = 0.71) with no-slip boundary conditions and observed 
periodic, quasiperiodic, and chaotic states for Rayleigh numbers between 6500 and 25000. Curry et al. (9f performed 
detailed DNS for P = 6.8 in 2D and observed oscillations with a single frequency at r ~ 50, and with two frequencies 
(quasiperiodic) at r ~ 290. They observed weak chaos beyond r ~ 290, and a periodic convective state after r ~ 800. 
Goldhirsch et al. [16j also simulated 2D RBC and observed complex behaviour. Zicnickc et al. [17| studied the effects 
of symmetries on 2D RBC using DNS, and observed periodic, quasiperiodic, and phase-locked states. 

Simulations of convection reveals various convective patterns, yet the large number of modes present in DNS 
obscures the origin of these patterns. Low-dimensional models, constructed using Galerkin projections, are very 
useful for investigating these issues. Some of the earlier investigations using low-dimensional models are by Lorenz [18| . 
Mclaughlin and Martin [19| and Curry [2fJ]. They observed various convective states including chaos through period 
doubling and quasiperiodic routes. Yahata [2l[ |22J constructed a 48-mode model that is inspired by the experimental 
configurations of Gollub and Benson [4| and reported quasiperiodic (Ruelle-Takens type) and period-doubling route 
to chaos. Recently Paul et al. [23| constructed a 30-mode model of 2D convection and studied it for P = 6.8 and 
P = 2\/2; they observed transition from periodic — > quasiperiodic — > phase-locked state — > chaos, similar to one of 
the experimental results of Gollub and Benson [J] . 

In this article we investigate the origin of various convective patterns using bifurcation diagrams constructed using 
DNS. We adopt pseudospectral method to simulate convection in water at room temperature (P = 6.8) in a box with 
aspect ratio of T x = 2V2. Bifurcation analysis requires a large number of computer runs. Since 3D simulations are 
computationally very expensive, we resort to 2D simulations. We obtain steady, periodic, quasiperiodic, and chaotic 
rolls as reported earlier by Curry et al. Q, Yahata [HI, and Mukutmoni and Yang [ll(. In addition, we observe that 
the flow becomes steady after chaos. Subsequently this steady state turns periodic and turbulent. We also observe 
coexistence of fixed points and chaotic attractors. We explain these features using bifurcation diagrams. 

The outline of the paper is as follows. In section II we describe the governing equations and the numerical method. 
Section III contains descriptions of various convective states resulting from various bifurcations. In section IV we 
discuss the dynamics of the chaotic state. We describe the properties of the large scale modes in section V. The last 
section contains conclusions. 

II. HYDRODYNAMIC SYSTEM AND NUMERICAL METHOD 

We consider a layer of Boussinesq fluid [l| of thickness d, kinematic viscosity v, thermal diffusivity k, and thermal 
expansion coefficient a confined between two stress-free and thermally conducting horizontal plates. An adverse 
temperature gradient ft = AT/d is imposed across the fluid layer, where AT is the temperature difference between 
the plates. The relevant hydrodynamic equations are nondimensionalized by choosing length scale as d, velocity scale 
as \J aftgd 2 , and temperature scale as AT, which yields 

i~~d 

(i) 

(2) 

(3) 

where v = (yi, 1)2,1)3) is the velocity fluctuation, 6 is the perturbations in the temperature field from the steady 
conduction state, R — agftd i /h i K is the Rayleigh number, P = v j k is the Prandtl number, g is the acceleration due to 
gravity, and z is the buoyancy direction. The two-dimensional rolls formed in the system are assumed to be parallel 
to the y axis. The top and bottom boundaries are considered to be stress free and perfectly conducting: 

^3 = d z V\ = d z v 2 = 8 = 0, at z = 0,1. (4) 

The velocity and temperature fields are periodic along the horizontal direction (x axis). 

The above set of equations ([T][3]) are solved numerically using a pseudospectral method [24j | in two dimensions under 
the above boundary conditions. We use Fourier basis functions for representation along the x direction, and sin or 




cos functions for representation along the z direction. The velocity component along the y-direction (i^) is zero. The 
expansion of the velocity and temperature fields are 
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(x,z,t) = \ 2U m Q n (t) exTp(imk c x)cos(mrz), 



v 2 (x,z,t) = 0, 

Vs(x, z,t) — y ^2W m Qn(t)exp(imk c x)sm(mrz), 

8(x,z,t) = y 29mOn(t)exp(imk c x)sm(nirz), (5) 

m,n 

where k c = ir/\/2. 

Various grid resolutions, 64 x 64, 128 x 128, 256 x 256, 512 x 512 have been used in the simulations. The aspect 
ratio of our simulations is T — 2\f2. Time stepping is carried out using the fourth-order Runge-Kutta (RK4) method 
with CFL scheme for choosing dt. The simulations are carried out till the system reaches a steady-state. 

We perform numerical simulations for 1.01 < r < 5 x 10 5 (R — 664 to 3.3 x 10 8 ). The thermal Prandtl number is 
chosen as P = 6.8 which is a typical value for water at room temperature. In the next section we will describe various 
convective states obtained in our simulations and investigate their origin using bifurcation diagrams. 

III. BIFURCATION DIAGRAMS FOR 2D RBC 

The number of variables of simulations on a N x N grid are TV 2 , which is quite large. It is impossible to understand 
the system dynamics in terms of all these variables. Fortunately numerical simulations reveal that some of the large- 
scale modes have large fraction of kinetic energy and entropy (J dx.\9\ 2 /2) up to the chaotic regime. Therefore we 
analyze various convective states of RBC using these large-scale modes, namely #101 and Wioi- In Figs. [T] and [2] we 
plot the absolute values of these modes as a function of r. For time-varying states like time-periodic and chaotic 
flows, we plot the minimum and maximum values of these variables. The above figures, also called the bifurcation 
diagrams, have been created from the results of around sixty DNS runs. In the following discussions we will describe 
the bifurcation diagrams in more detail. 

The phase space or state space of a system provides valuable information on the dynamics of that system. The 
phase space of our dynamical system is N 2 dimensional, which is impossible to visualize. In the following discussion 
we will show only a projection of the phase space onto a subspace made by some of the energetic modes. In Fig. [3] we 
plot a phase space projection onto Q(W±oi) — ^(#101) plane, where 9 stands for the imaginary part of the argument. 

The RBC system is in a conduction state for r < 1. At r = 1, the conductive state becomes unstable, and a 
convective state is born via a supercritical pitchfork bifurcation. A new stable solution, the stationary straight rolls, 
represented by blue diamonds in Figs. Q] and [5] emerge from the bifurcation point. The stationary straight rolls exist till 
r ~ 80. A phase space projection of the fixed points is shown in Fig. G3a). Various secondary bifurcations take place 
after this primary bifurcation. Before we proceed to discuss various convective states, we make an important remark 
regarding the fixed point solutions and related states. When the initial conditions of the velocity and temperature 
Fourier modes are chosen as imaginary (except only for #oon which is purely real) , the fixed points are purely imaginary. 
However, when we choose complex Fourier modes as initial conditions, we get fixed points with complex Fourier modes, 
but the absolute values of all the Fourier modes match with the magnitude of the above mentioned imaginary values. 
This feature is a consequence of the translational invariance of the solutions in the x-direction due to the periodicity 
of the box along x. This manifests itself as a phase shift of the various Fourier modes. Therefore a purely imaginary 
mode can become complex with the same magnitude. This symmetry plays an important role in the chaotic regime 
as well; these results will be described in our later discussions. 

Near r ~ 80, the system bifurcates from a steady convection state to a periodic state (a limit cycle in the phase 
space) through a Hopf bifurcation. The extremum values of |Wioi| and |#ioi| are shown in Figs. [T] and [2] respectively 
using red dots. A sketch of the unstable fixed point solutions obtained using extrapolation is also shown in the 
bifurcation diagrams as a dashed blue curve. Figure [3Jb) exhibits a phase space projection of a limit cycle on the 
3(^101) - 3(#ioi) plane for r = 100. The power spectrum of the time series of the Fourier mode 5(Wioi) for the same 
r is shown in Fig.[4ja) clearly depicting sharp peaks at / = f\ (dominant) and 2/i (superharmonic) with /1 sw 0.2624. 
As r is increased further, the periodic state undergoes a period-doubling bifurcation at around r ~ 350, similar to 
that observed by Mukutmoni and Yang (ll) . As a consequence of this bifurcation, the phase space projection of the 
limit cycle on the Q'(Wioi) - 3(6*ioi) plane shows a second loop as depicted in Fig.[3][c) for r = 350. For this state, the 
power spectrum acquires a peak at /i/2 along with the dominant peak at f\ as shown in Fig. [4jb). Superharmonic 



frequencies are also present in Fig. @|b) due to the nonlinearity in the system. We note at this point that the period-2 
behavior is not very apparent in the bifurcation diagrams (Figs. Q]and[2]) f° r the absolute values of the modes W^oi 
and ^loi- This is possibly due to some cancellations while taking the absolute values of the modes. 

The period-2 limit cycles described above become unstable at around r ~ 500, and another independent frequency is 
born through a Neimark-Sacker bifurcation. The power spectrum of 5(Wicu), shown in Fig.|4jc) for r = 500, indicates 
three dominant peaks at frequencies /i/2 ~ 0.158, fa ~ 0.316 and fa ~ 0.0075. The first two frequencies correspond 
to the period-2 limit cycle, while the third frequency is an irrational multiple of the first two (fi/fa ~ 42.2362). As 
a consequence, the RBC state for r = 500 is quasiperiodic, and the phase space trajectories fill a torus in higher 
dimensions, whose projection on a plane is shown in Fig. [3jd). The quasiperiodic flows, observed for 500 < r < 725, 
are represented by green dots in the bifurcation diagrams (Figs. [TJ [5]). 

At r ~ 750, the quasiperiodic state eventually bifurcates to a chaotic state that persists till r ~ 849. In Figs. [TJ and [2] 
the chaotic flows are illustrated as pink dots. A phase space projection of a chaotic state at r = 750 is shown in Fig. (3Je) . 
The power spectrum for this state, shown in (Fig. HJd)), shows a broadband power spectra, a characteristic of a chaotic 
solution. This chaotic attractor (Fig. EJe)) observed in the range 750 < r < 840, is born through a quasiperiodic 
route to chaos. This transition could either be through phase- locking or through Ruelle-Takens scenario ( |25l. 1 26jj ) . 
which can be ascertained only by performing simulations for many r values in the neighbourhood. DNS cannot help 
us ascertain which of the two scenario is applicable since simulations for many r values in the neighbourhood are not 
possible. Yet close similarity of DNS results with the corresponding low-dimensional calculations [23J suggests that 
the transition is possibly through the phase-locking route. 

In the band 750 < r < 840 we observe two unconnected chaotic attractors that lie in the first and the third 
quadrants of the phase space projection. That is, both Qf(Wioi) and 9f(0ioi) either take positive values or negative 
values. At later r values, for 841 < r < 849, these two disjoint attractors merge through an "attractor merging 
crisis" [2g, [27| to create a larger chaotic attractor shown in Fig. |3Jf ) . However the increase in the size of the chaotic 
attractor is not apparent in the bifurcation diagrams (Figs. [1] and [2]) due to the use of absolute values of these modes 
for depiction. 

As r is increased further, we observe two coexisting attractors in the window of 846 < r < 849. For some initial 
conditions we observe chaos, and for some others we find stable fixed points. For r = 846, the coexisting attractors 
in a phase space projection are shown in Fig. [5] wherein the dots represent stable fixed points, while the trajectory 
corresponds to a chaotic solution. This feature can be understood as follows. The unstable fixed points (dashed blue 
line in the bifurcation diagram) undergoes an inverse subcritical Hopf bifurcation at r ~ 845. As a result, an unstable 
limit cycle is born (dashed red lines in the inset of Fig. [TJ), which coexists with the stable fixed points (the solid blue 
line in the inset) for r > 845. As a result of this bifurcation, stable fixed points and a chaotic attractor are observed 
to coexist for 846 < r < 849 with the unstable limit cycle (along with its stable manifold) approximately forming the 
basin boundary between the two. 

At r ~ 849, the chaotic attractor disappears, and we observe only stable fixed points as depicted in Fig. [3Jg) for 
r = 1500. The disappearance of the chaotic attractor is due to a "boundary crisis" [26|,[28| where the chaotic attractor 
collides with the unstable limit cycle. Figure |5] shows the time series of K(#ioi) and 3(#ioi) for r = 847 and 850. The 
system is chaotic for r = 847, however for r = 850 it settles down to a stable fixed point after chaotic transients. A 
close similarity between the long transients for r = 850 and the chaotic time series for r — 847 provides an evidence 
for the occurrence of a "boundary crisis" . 

The new fixed points continue to remain stable for a significantly large range of r values. At higher r values, the 
fixed points bifurcate to limit cycles through a Hopf bifurcation. A phase space projection of a limit cycle for r = 7000 
is shown in Fig. [3jh). At later r values, the limit cycles bifurcate to a turbulent state, similar to that reported by 
Vincent and Yuen [20|. A sample phase space projection for this range is shown in Fig. \%li) for r = 5 x 10 5 . The 
properties of turbulent state has been explored earlier by Vincent and Yuen [29| for 10 8 < r < 10 14 , and they are not 
discussed in this paper. Our bifurcation diagrams shown in Figs. [1] and [5] are quite detailed till r = 1600. 

The quasiperiodic route to chaos for 2D RBC presented here is in genera l ag reement with those described earlier 
in some experiments [J], simulations ISH121 fl5j ■ and low-dimensional models 23]. There is a similarity with the RBC 
simulations of Mukutmoni and Yang |ll| . who observed a period-two state but subsequent route to chaos through 
quasiperiodicity We report for the first time stable fixed points beyond chaos, and coexistent stable fixed points and 
chaotic states. We also observe several "crisis" for the first time in 2D RBC simulations: first, smaller attractors 
corresponding to Fig. \5j[e) merge to form a larger attractor shown in Fig. (3]Jf ) through an "attractor-merging crisis" , 
and second, the disappearance of the chaotic attractor through a "boundary crisis" . 

In the present paper we focus on explaining the origin of various convective patterns using detailed bifurcation 



analysis. Recently Paul et at 23j performed similar studies using a 30-mode model that was derived using a Galerkin 
projection of 2D RBC on appropriate large-scale modes. The features presented here are quite similar to the bifurcation 
scenario presented for the model (Figs. [7] and H]). There are however some differences, e.g., the low-dimensional model 
has limit cycles and a chaotic attractor as coexisting attractors, while the DNS shows fixed points and a chaotic 



attractor as coexisting attractors. Another difference stems from the fact that for a given pattern, the r values for 
2D DNS is larger than the corresponding model values. 

Experiments, DNS, and low-dimensional models reveal strong dependence of the convective patterns and emergence 
of chaos on the Prandtl number and the aspect ratio. Our DNS results are for P = 6.8 and aspect ratio of 2\/2- 
We need to perform more extensive DNS to come up with a comprehensive picture of large-P convection. After the 
discussions on bifurcation diagrams, we describe the dynamics of chaotic states in the next section. 

IV. DYNAMICS OF THE CHAOTIC STATES 

As described in the previous section, chaos appears in 2D RBC for P = 6.8 and V = 2\/2 through a quasiperiodic 
route. The chaotic attractor for 750 < r < 840 is confined in a narrow region of the phase space, while the region of 
the chaotic attractor for 841 < r < 849 is much larger. This difference between the sizes of the chaotic attractors can 
be understood by projecting the phase space on ^(W^ioi) — ^(Wioi) plane (see Fig. [9]). The figure in the left panel 
corresponds to the narrower chaotic attractor for r = 790, similar to Fig. G3e). In the left panel, the phase space 
trajectories for r = 790 wander about chaotically in the region K(Wioi) = (-0.019 : 0.0024), ^(Wioi) = (-0.10 : 
—0.077), with the phase of the Fourier mode Wioi confined in a narrow band of (1.437T : 1.5l7r) radians. Note that 
the phase of the Fourier mode W\oi could take any mean value due to azimuthal symmetry discussed in the previous 
section. 

For 841 < r < 849, the chaotic attractor increases in size as a result of an "attractor-merging crisis" . A sample of 
such a chaotic attractor is shown in the right panel of Fig. [9] for r — 841. The system evolves in such a way that the 
phase of the mode W\q\ varies in a narrow band for a while, and then the system jumps abruptly to another band of 
phases. Here, the phase of the mode W\q\ can take any value between and 2ir unlike the smaller chaotic attractor 
described above. The region occupied by the phase space attractor at r = 790 is confined within the dashed box of 
the right panel of Fig. [9] 

The above mentioned abrupt changes in the phase of the Fourier mode Wioi correspond to sudden or chaotic 
movements of the convective rolls along the horizontal direction [30] ■ We show this feature in Fig. QJJ where we have 
drawn the snapshots of the flow patterns for r = 841. The three frames of Fig. [TU] represent t — 14000, just before 
the rolls start moving, t = 15000, during the movement of the roll, and t — 16000, when the roll has moved. As 
mentioned earlier, the movements of the rolls are quite abrupt. A movie depicting the movement of the roll during 
this time interval is available at [311 ] . 

The travelling rolls are also related to the reorientations of convective structures and wind reversals observed in 
experiments and simulations 32-34]. We also observe interesting patterns in the energy content of the dominant 
Fourier modes during these abrupt events. The time series of |Wioi| shown in Fig. [TT] shows a decrease in the values 
of iPV'ioil up to around 20% of its maximum value when the phase of the Fourier mode changes abruptly. The value 
of IW202I however increases suddenly during these events. These features of the energetics of the Fourier modes are 
similar to those observed in the DNS of wind reversals in a cylindrical container [34J ] . 

As described earlier, convective states like steady, periodic, quasiperiodic, and chaotic rolls can be conveniently 
described using the large scale or small wavenumber Fourier modes. Variation of the magnitudes of these Fourier 
modes with r follows some interesting patterns that will be described in the next section. 

V. PROPERTIES OF THE LARGE-SCALE MODES 

In Fig.[12]we plot the values |Wioi| and |(9ioi| for the stable fixed points as a function of (r — 1). We observe that 
|Wioi I ~ ( r — l) ' 62 for r up to 1000. However |#ioi| shows two distinct scaling: ~ (r — l) 027 for 1 < r < 2, and 
~ (r — 1)~ - 34 for 2 < r < 1000. We also compare the values of |Wioi| and |0ioi| computed using DNS, the 30-mode 
model of Paul et al. [23|, and the Lorenz model [26[ for P — 6.8. The entries listed in Table Q] indicate that the results 
of the 30-mode model matches quite well with DNS for r < 10, while Lorenz model works well for r only up to 4. 
For r far away from the onset, many Fourier modes get excited, and the low-dimensional models are not expected to 
match with DNS quantitatively. However the energetic Fourier modes provide fairly good qualitative description as 
visible from the significant similarities between the DNS bifurcation diagrams (Figs. [T] and [2]) and the corresponding 
diagrams obtained from the low-dimensional model (Figs. [7] and [5]). 

Another quantity of interest related to the large-scale convection is the Nusselt number, which is the ratio of the 
total heat flux and the conductive heat flux. Experiments, numerical simulations, and phenomenological theories 
project a power law behaviour of the Nusselt number as a function of the Rayleigh number [35j. The predicted 
spectral indices for various conditions are 1/4, 2/7, 3/10, 1/3, 1/2, etc. (35(. The Nusselt numbers Nu for our various 
DNS runs are shown in Fig. [TBI where we plot Nu as a function of (r — 1). For r > 2, Nu ~ (r — i)°- 33 ± 001 with 
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TABLE I: Comparison of values of |VFioi| and |#ioi| from the DNS, the 30-mode model and the Lorenz model. The values 
match closely for low r. With higher r the values start to differ significantly. 

an observable deviation from the fit for 100 < r < 600. These results are in good agreement with the experimental 
findings of Cioni [32| and Niemela et al. [36| . as well as many numerical simulations, e.g., by Moore and Weiss |L4J |. 
This agreement also acts as a validation of our numerical code. 



VI. CONCLUSIONS 



In this paper we perform direct numerical simulations of two-dimensional RBC for P = 6.8, aspect ratio T = 2\/2, 
and the reduced Rayleigh number r up to 5 x 10 5 (R = 3.3 x 10 7 ). Simulations reveal various convective states like 
steady, periodic, quasiperiodic, and chaotic rolls. Using simulation results we construct bifurcation diagrams that help 
us understand the origin of the observed convective patterns. The steady convective rolls are born through a pitchfork 
bifurcation at r — 1. At around r ~ 80, time-periodic rolls appear through a Hopf bifurcation. The periodic state 
bifurcates to period-2 state through a period-doubling bifurcation. Through a Niemark-S acker bifurcation, these 
period-2 rolls turn into quasiperiodic rolls which in turn becomes chaotic through a quasiperiodic route to chaos. 
The chaotic attractors undergo an "attractor-merging crisis" to generate a larger chaotic attractor involving abrupt 
changes in the phase of the Fourier modes leading to abrupt motion of the convective rolls. In our DNS we also 
observe coexistence of stable fixed points and a chaotic attractor as a result of an inverse subcritical Hopf bifurcation. 
After a while, only stable fixed points survive due to the disappearance of the chaotic attractor through a "boundary 
crisis" . These fixed points subsequently bifurcate to periodic and chaotic states through successive bifurcations. 

Several of the above mentioned convective structures have been observed earlier in some experiments [J, [5J, |8j, 
simulations |9l-[r2L Ha ], and low-dimensional models [2l], (23|. However, we observed for the first time stable fixed 
points beyond chaos, coexistent stable fixed points and chaotic state, and several occurrences of "crisis" . The above 
mentioned bifurcation diagrams have similarities with those derived for Paul et o/.'s 30-mode model for the same 
Prandtl number [23j. 

Experiments of Gollub and Benson [4] exhibit other kinds of patterns and chaos for different sets of Prandtl numbers 
(2.5, 5) and aspect ratios (2.4, 3.5). We are exploring these ranges of Prandtl numbers and aspect ratios through 
both DNS and low-dimensional models and hope to achieve a more comprehensive picture of large-P convection in 
future. The results in our paper illustrate the usefulness of bifurcation diagrams in understanding the origin of various 
convective patterns. 
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FIG. 1: Bifurcation diagram |#ioi| vs. r computed using DNS results. The system undergoes a pitchfork bifurcation at r = 1. 
The new stable fixed points (blue diamonds) undergo a supercritical Hopf bifurcation at r ~ 80 yielding stable time-periodic 
flows (red circles). Quasiperiodic solutions are born at r ~ 500. Subsequently chaos is observed in the band 750 < r < 849. 
Coexisting attractors, stable fixed points and a chaotic attractor, are observed for 846 < r < 849 as shown in the right inset. 
Subsequently the chaotic attractor disappears through a "boundary crisis" at r = 850, and only stable fixed points survive. 
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FIG. 2: Bifurcation diagram |Wioi| vs. r computed using DNS results. The symbols used for various states and their 
interpretations are same as Fig. [T] 
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FIG. 3: For P = 6.8, phase space projections on SJ(Wioi) - Ss(6*ioi) plane of various convective states: steady convective rolls 
for r = 70 (a); time-periodic state for r — 100 (b); time-periodic roll with period-2 behavior for r = 350 (c); quasiperiodic state 
for r — 500 (d); chaotic states for r = 750 (e) and for r = 841 (f); reemergence of steady convective rolls for r — 1500 (g); 
periodic state for r — 7000 (h); turbulent state for r = 5 x 10 J (i). 
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FIG. 4: Power spectrum of Sy(Wioi) for various convective states: (a) for periodic state with dominant peak at /i = 0.2632 
(r = 100); (b) for period-2 state with dominant peaks at /i/2 and /i (r = 350); (c) for quasiperiodic state with dominant 
peaks at /i/2, /i, and fa with /1//2 ~ 42.2362 (r = 500); (d) for chaotic state at r — 750 with broadband spectrum. 
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FIG. 5: For r = 846, a phase space projection of coexisting attractors, stable fixed points (red dots) and a chaotic attractor 
(blue trajectory). 
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FIG. 6: Time series of 5R(#i i) and Q^ioi) for r — 847 and r — 850. At r = 847 the system is chaotic, while at r = 850 the 
system settles down to a fixed point solution after chaotic transients. 
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FIG. 7: Bifurcation diagram computed using the 30-mode model of Paul et al. [23]. FP (blue curve) is the steady convection 
roll, 'OS' (red curve) is the time-periodic roll, 'QP' (green patch) is the quasiperiodic roll, and 'CH' (pink patch) is the chaotic 
state. 'NS' indicates the Neimark-Sacker bifurcation point. A window of periodic and quasiperiodic states is observed in the 
band of r = 46.2 : 48.4. This bifurcation diagram has certain similarities with that computed using DNS results (Fig. [2]). 
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FIG. 8: Bifurcation diagram for |#ioi| computed using the 30-mode model of Paul et al. [23]. Various bifurcations are similar 
to as in FIG. [7] However this bifurcation diagram is for a different initial condition, which does not show a window of periodic 
and quasiperiodic states after chaos. This illustrates the coexistence of periodic and chaotic states in the range r = 46.2 : 47.4. 
For comparison with DNS bifurcation diagram, see Fig. [1} . 
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FIG. 9: Phase space projection of the chaotic attractor on the 5ft(Wioi) — Ss(Wioi) plane for (a) r — 790 for which the phase 
of the mode Wioi lies in the band of (1.437T : 1.5l7r); (b) r — 841 for which the phase of VKioi takes all values in (0, 2-k). The 
range of the chaotic attractor of figure (a) is depicted by the dashed box of figure (b). 
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FIG. 10: The snapshots of the chaotic travelling rolls at r = 841 and P = 6.8. The convective rolls move to the left. 
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FIG. 11: Time series of [Wioi| and IW202I modes for r — 841. | W101 1 value decreases sharply, while IW202I rises when the phase 
of W101 changes abruptly. 
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FIG. 12: Variation of |Wioi| and |0ioi| as a function of (r — 1) obtained from the DNS results. |Wioi| ~ (r — l) 0,62 , while 
\8wi\ ~ (r- l) 027 for r < 2, and |0ioi| ~ (r - l)~ - 34 for 2 < r < 1000. 
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FIG. 13: Variation of Nusselt number Nu as a function of (r-1). Nu ~ (r - i)°-«± - 01 for r > 2 



